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We calculate for the first time the complex potential between a heavy quark and antiquark at finite 
temperature across the deconfinement transition in lattice QCD. The real and imaginary part of the 
potential at each separation distance r is obtained from the spectral function of the thermal Wilson 
loop. We confirm the existence of an imaginary part above the critical temperature Tc, which grows 
as a function of r and underscores the importance of collisions with the gluonic environment for the 
melting of heavy quarkonia in the quark-gluon-plasma. 



Heavy-quark bound states (QQ) are essential tools 
in the experimental and theoretical investigation of the 
high temperature state of QCD matter, the quark-gluon 
plasma (QGP) Q. In particular, the suppression of 
heavy quarkonia such as J 'ftp (the ground state of cc in 
the vector channel) in relativistic heavy-ion collision ex- 
periments at the Relativistic Heavy Ion Collider Q and 
at the Large Hadron Collider [3J, provides us with an 
intriguing signal of the QGP. 

In order to extract the physics of the QGP from heavy 
quarkonia, an intuitive as well as quantitative under- 
standing of the involved physics is necessary: In the 
classic approach of using a Schrodinger equation with 
a phenomenological QQ potential, Debye screening of 
color charges is responsible for the melting of QQ bound 
states above the deconfinement temperature Tq [J]. A 
perturbative evaluation of the thermal Wilson loop us- 
ing the hard-thermal-loop resummation technique shows 
that the static potential at high temperature has an 
imaginary part induced by Landau damping Q . An 
approach based on heavy-quark effective field theory at 
finite temperature indicates an additional contribution to 
the imaginary part of the static potential [7j . In addition, 
lattice QCD simulations of the spectral function of char- 
monia at finite temperature show that J/ip may survive 
even up to 1.5Tc However, the connection among 
the different approaches is not yet clearly understood. 

The main purpose of this Letter is to unify the ap- 
proach based on the nonrelativistic Schrodinger equa- 
tion and that based on a spectral decomposition of the 
QQ correlator. Utilizing lattice QCD simulations of the 
medium surrounding the QQ, we extract a nonperturba- 
tive potential at any temperature, especially in the phe- 
nomenologically important and hitherto inaccessible re- 
gion around Tc- 

The starting point for deriving the effective 
Schrodinger equation at finite T is the real-time 
forward QQ corr elator (t > 0) [BJ 



D 



•(r,t) = (M(x,y,t)Mt(x,y,0) 



where M(x } y, t) — Q(x,t)jU[x,y]Q(y,t) with 7 denot- 
ing gamma matrices (7 = 7^ for J/ip) and L7[x, y] being 
a straight Wilson line connecting (x, i) and (y,t). The 
relative coordinate is defined as r = x — y. The angu- 
lar brackets impliy the thermal average. In the nonrel- 
ativistic limit, where the heavy-quark mass rriQ is much 
larger than any other scale of the system, one may rewrite 
Eq. ([IJ in terms of a path integral over the positions and 
momenta of the heavy quarks: 
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with the Wilson loop, Wr = ^exp 

Here T denotes a closed loop in Minkowski space-time, 
defined by the path of the quark and antiquark. This 
formula is a finite temperature generalization of the result 
in Ref. [9J and is valid up to 0(\/mQ). 

From here on, we focus on the leading part of the po- 
tential, obtained in the heavy-quark limit toq — » 00. In 
this case, it is enough to consider a rectangular path for 
the Wilson loop Wr — > Wo ( r = l z i — z 2 1 , t) with a Fourier 
decomposition; 
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Equation (fT]) for ttiq — > 00 also has a spectral decompo- 
sition composed of three parts, p\ oop , p s t ap ic and phandi c : 
one can show that pn(r, uS) in Eq.Q corresponds to pi 00 p 
and is positive semidefinite for all r and to [10(. 

Through Eq. (J3J) we are led to the following nonpertur- 
bative definition of the in-medium potential 



id t W a (r,t) _ f duue-^pufau) 



Vp(r,t). (4) 



Indeed, the solution of Eq. Q yields the transfer matrix 



(1) compatible relation Wp(r, t) = exp — i J Q Vo(r, s)ds 
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FIG. 1. (Left) Lattice data for the Euclidean thermal Wilson-loop as a function of the imaginary time r at T/Tc = 1.17. 
The solid lines show the results obtained after the MEM reconstruction of the spectral function. (Right) Spectral functions 
and their prior dependence obtained by the MEM with the thermal Wilson- loop data at r/a a = 1,4. Positive cj regions and 
negative ui regions are plotted separately with a different vertical scale. 



and thus leads to an in-medium Schrodinger equation 
id t D> R (r,t)= (-%- + V n (r,tj)D> R (r,t). (5) 



In the following we study the peaks of the spectral 
function in Eq.Q around uj = 0, since it is the physics 
at late times t that determines the potential. Based on 
Eq.dU, the real part of the potential is thus given by 
the peak position, whereas the imaginary part is related 
to the peak width: For example, if p\j(r, uj) exhibits a 
Breit-Wigner peak with position u>o(r) and width Tq{t), 
we obtain V D (r,t) -> V^ BW >(r) = ui (r) + iT (r), while 
if it contains a single Gaussian peak, we have V\j (r, t) — > 
VW{r,t) = wo(r) + ir§(r)t. 

A nonperturbative method to determine /Qrj(r, uj) and 
hence Vp(r, t) for all T is based on lattice QCD simula- 
tions. We first make an analytic continuation of Eq.Q 
to imaginary time and connect the Euclidean thermal 
Wilson-loop to the spectral function through the Laplace 
transform, W^(r 7 r) = J duie~^ T p<^{r, uj). This quan- 
tity is purely real and not symmetric under the reflection 
t /3 — t by definition. As shown in Ref. [lOj . pn(r, uj) 
in a finite volume allows for a spectral decomposition in 
terms of delta-function peaks: 



, P n (T)\A 1hn ,(r)\ 2 5(uj-E% Q (r)+E n ). (6) 



Here E n is the energy eigenvalue of a state \n) without 
QQ, while E®®(r) is the eigenvalue (relative to 2mg) of 
a state \nf; r) including a QQ separated by a distance r. 
The matrix element and the Boltzmann factor are defined 
as A nn ,(r) = £ y (n|M(y + r,y,t = 0)|n') and P n (T) = 
e -E„/T j [J2 n e~ En / T ] , respectively. The peak positions 

of Eq.((6|) do not depend on T (since E n and E®,® (r) are 
T- independent), while the pole residues P n {T)\A n , n i (r)\ 2 
depend on T due to the Boltzmann factor. The various 



possible combinations of n' and n in pu(r, u) at finite T 
lead to a number of peaks as large as N x N, if pp(r, uj) 
at T — has N peaks. 

The position and width of a physical resonance in an 
infinite volume are obtained from the envelope of a bunch 
of delta- function peaks in Eq. © . A suitable technique to 
extract such envelopes from lattice QCD data of Wpj (r, r) 
is the maximum entropy method (MEM) [ll|. In this 
Letter, we utilize a high precision MEM with an extended 
search space developed by one of the present authors [l2| ■ 
By separating the size of the search space from the num- 
ber of data points, it allows a consistent comparison of 
spectra even if the temporal extent of the underlying lat- 
tice is changed. 

The stability of the position and width of the reso- 
nance reconstructed by the MEM can be checked through 
changing the amplitude of the prior and by increasing the 
number of temporal data points ll| . In addition, the ef- 
fect of a discrete lattice spectrum on the position and 
width of the resonance can be checked by controlling the 
number of peaks and their relative separation under the 
variation of the lattice spacing a and the lattice size L 
(see e.g. [HI for the variation of a). 

We perform quenched lattice QCD simulations of 
W®(r 7 T) by using the simple plaquette gauge action 
on an anisotropic 20 3 x N T lattice. The anisotropy ra- 
tio between the the spatial and temporal lattice spac- 
ing is taken to be £ = a a /a T — 4. We fix a a to 
be 0.097fm ((3 = 6.1) and adopt the temporal lattice 
sizes N T = 36, 24, 12 which correspond to temperatures 
T/T c = 0.78, 1.17, 2.33 with T c ^ 290MeV[2|. Our spa- 
tial lattice size L ~ 2fm can accommodate the character- 
istic J/i/j scale of rj/^p ~ 0.5fm. After collecting the lat- 
tice data, we carry out the MEM over a frequency interval 
lui = [Umm ! w max ] of N u = 1500 points. Our choice cor- 
responds to I u ~ [-21GeV,42GeV] at N T = 24. We use 
a prior distribution of the form m(ui) = , motivated 
by the canonical dimension of pair, uj). The parameter 
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FIG. 2. (Left) The real part of the in-medium heavy-quark potential for three different temperatures across Tc with a a = 0.0971 
fm (solid squares, circles and triangles). Solid and dashed lines are the color-singlet free energies F x (r) in the Coulomb gauge 
at corresponding temperatures. The result for the finer lattice (a CT = 0.039 fm at T /Tc = 2.33 is also shown for comparison 
(solid diamonds). (Right) The imaginary part of the potential, obtained from the width of the lowest lying peak of the spectral 
function. The same symbols are used as the left figure. 



uiq(> 0) is fixed by setting the amplitude mo = m(w m i„) 
at the smallest frequency in I u . 

In FigQJleft) we pl°t typical Wilson-loop data as a 
function of r for various distances r at T = l.YITc- Spec- 
tral functions obtained by the MEM at r — a a ~ 0.1 fm 
and r — 4a CT ~ 0.39 fm are plotted in FigQJright). The 
falloff of W^(r,r) for small and intermediate values of 
t in the left figure corresponds to the peaks located in 
the u> > region seen in the right figure: they arise from 
gluonic interactions between the temporal Wilson lines. 
The upward trend of W§ (r, r) around r = /3 on the other 
hand is induced by extremely small structures located in 
the lo < region: These arise from short distance gluon 
interactions connecting the spatial Wilson lines across the 
compactified temporal axis. Although the spectral func- 
tion in the negative uj region is important for reproducing 



/3) due to its exponential "enhancement" by 



the Laplace transform, its effect on V\j(r,t) in Eq.(j4]) is 
negligible due to its extremely small residue. 

The real and imaginary parts of the potential are ob- 
tained from fitting the lowest lying peak with a Brcit- 
Wigner and Gaussian shape. The identical results for 
position and width are shown in Fig|2]deft) and (right), 
respectively. The error bars are obtained from the vari- 
ance in the peak structure between different choices for 
the amplitude of the prior distribution. The interval of 
10~ 2 > mo > 10~ 6 is chosen to span as many orders 
without introducing numerical instabilities in the min- 
imization process. We have checked that the standard 
MEM error, estimated from the stability of the spectral 
function given mo ll|, is much smaller than the error 
from the variation of mo . 

T = 0.78T C : Rc[Vb(r)] denoted by the filled squares 
in FigHJleft) is found to show a linearly rising potential 
at long distances, a result consistent with quark confine- 
ment below Tq- To study the short distance behavior 
of the potential, we have also carried out simulations 



with a 2.5 times smaller lattice spacing ((3 = 7, £ = 
4, a a = 0.039fm, 20 3 x 96) at the same T/T c ; the result- 
ing Re[Vn(r)j can be fitted well by a Coulomb 4- linear 
potential [151 ]. Furthermore, our results agree with the 
color-singlet free energies F 1 (r) in the Coulomb gauge 
(solid black line in the left panel) at all measured dis- 
tances r. Note that Tn(r) appears to be small and is 
consistent with zero within the statistical and systematic 
errors; see the filled squares in Fig|2]Jright). 

T = 1.12T C : At this temperature Re[Vb(r)] has ap- 
parently frozen at around the same strength found at 
T = 0.78Tc, as shown by the filled circles in the left 
panel. This is in contrast to the behavior of F 1 (r) at the 
same T (the long dashed curve in the left panel) which 
exhibits a significant thermal screening. As for rrj(r), 
there is a tendency to develop a nonzero value, which 
grows as r increases; see the filled circles in the right 
panel. 

T = 2.33Tc: At the highest available temperature, 
Re[Vb(r)] and T\j(r) exhibit a strong rise as a func- 
tion of r, shown by the filled triangles. To test the 
effects of a small number of temporal data points N~ 
[IH and also whether the discrete lattice spectrum [l3| 
might blur the MEM image and thus might lead to an 
artificial broad peak, we compare the results of (j3 = 
6.1, a a = 0.097fm,20 3 x 12) and those of (0 = 7,f = 
4, d(j = 0.039fm, 20 3 x 32): The latter has more temporal 
data points; hence, it has finer (coarser) resolution at low 
(high) frequency. The results shown by solid diamonds 
in Figj2{left) and (right) are consistent with the solid tri- 
angles. This crosscheck indicates that the width broad- 
ening observed from the MEM is a physical effect [l6j |. 
We need, however, further systematic studies to confirm 
this point by changing N T , a and L independently. 

Although we construct our potential by starting from 
the mesonic operator M(x, y, t) with a straight Wilson 
line in Eq.([T]), one has the freedom to choose the operator 
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FIG. 3. The real part (upper panel) and the imaginary part 
(lower panel) of the potential obtained from the Wilson lines 

M differently and, hence, the forward correlator D > and 
the potential V. On the other hand, observable quanti- 
ties, such as the dilepton emission rate must be indepen- 
dent of such differences 171. Therefore, there must be a 



trade-off between the real and the imaginary part of the 
potential to leave the observables unchanged. To study 
this point, we consider an operator M with U[x, y] = 1 in 
the Coulomb-gauge. In this case, we have "Wilson lines" 
without the spatial link W\\{r,t). In Fig 02 we show the 
real and imaginary parts of the potential (Vj|) obtained 
from W®(r,r). In the confinement phase below Tc, V\\ 
and Vq agree quite well. On the other hand, in the de- 
confinement phase, both Re[Vji(r)] and T\\{r) exhibit a 
less pronounced rise in r. We observe that a weaker real 
part is accompanied by a weaker spectral width, which 
could be a sign of the trade-off mentioned above. This 
mechanism has to be made quantitative in future studies 
by solving the time dependent Schrddinger equation for 
an initial QQ wave packet entering the QGP. 

We have presented a nonperturbative derivation of the 
Schrodinger equation for heavy quarkonia and a first eval- 
uation of the corresponding complex in-medium poten- 
tial, based on quenched lattice QCD. Our numerical re- 
sults show that, even though the potential agrees with 
the color-singlet free energies below the phase transition, 
the correct physics above Tc can be obtained only if the 
real and imaginary part are taken into account together. 
The temperature insensitive real part around Tc suggests 
furthermore that the growth of the imaginary part, i.e. 
an increasing number of collisions with the medium, may 
play a more important role to destabilize QQ than the 
screening effects [l8j]. We have also discussed a possible 
mechanism to obtain the relevant physics independent of 
a particular choice of the underlying operator by balanc- 
ing the real and imaginary part. Our complex poten- 
tial opens up new possibilities to study the dynamics of 
the QGP transition by providing first-principles input to 
nonrelativistic real-time simulations, going beyond both 



models and perturbation theory. 

Our ongoing work aims at full QCD simulations with 
dynamical fermions, since these additional degrees of 
freedom may affect both the real and imaginary part of 
the complex potential substantially. In addition, larger 
and finer lattices are needed in order to assess the relative 
significance of Debye screening vs the collisional effects 
from short distance to long distance in more detail. 
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